Autophagy-related biomarkers in preeclampsia: the underlying mechanism, correlation to the immune microenvironment and drug screening

Background Preeclampsia is a life-threatening disease of pregnancy that lacks effective pharmaceuticals which can target its pathogenesis. Since preeclampsia involves complex pathological processes, including autophagy, this study aims to explore autophagy-related mechanisms of preeclampsia and to screen potential drugs. Methods Firstly, the datasets GSE75010, GSE24129, GSE66273, and autophagic genes lists were downloaded from public databases. Then, a weighted gene co-expression network analysis (WGCNA) was applied to filter autophagic-related hub genes of preeclampsia. The differential expression levels of the hub genes were validated with datasets GSE24129 and GSE66273. Next, the GO and KEGG enrichment, protein-protein interacting (PPI) network, as well as the downstream pathways was analyzed via the starBase, STRING and Cytoscape to determine the functions and regulatory network of the hub genes. Additionally, the immune microenvironment of preeclampsia was investigated by the CIBERSORTX database. Finally, three herb ingredients, berberine, baicalein, and luteolin were screened by molecular docking in comparison to pravastatin, metformin, and aspirin, to predict potential drugs for treating preeclampsia. Results A total of 54 autophagy-related genes were filtered by WGCNA. After filtering with |GS| > 0.5 and |MM| > 0.8, three hub genes, namely PKM, LEP, and HK2, were identified and validated. Among these genes, PKM and LEP were overexpressed in women older than 35 years old ( p<0.05; p<0.05); the expression of PKM, LEP, and HK2 differed remarkably in women with different BMI (all p<0.05); PKM overexpressed in women with hypertension (p<0.05). The regulatory network of hub genes demonstrated that they were mainly enriched in metabolic pathways, including the AMPK signaling pathway, glucagon signaling pathway, adipocytokine signaling pathway, and central carbon metabolism. Then, immune microenvironment analysis turned out that M2 macrophages were reduced in preeclampsia women (p<0.0001) and were negatively correlated with the expression of PKM (r=-0.2, p<0.05), LEP (r=-0.4, p<0.0001), and HK2 (r=-0.3, p<0.001). Lastly, molecular docking showed baicalein and luteolin could bind intimately to hub genes. Conclusion PKM, LEP, and HK2 could be promising biomarkers for preeclampsia, which might regulate the pathogenesis of preeclampsia via metabolism pathways and immune microenvironment. Baicalein and luteolin could be potential therapeutics for preeclampsia. Supplementary Information The online version contains supplementary material available at 10.1186/s12884-023-06211-2.


Introduction
Preeclampsia (PE) is a multi-system disorder that classically presents with newly onset hypertension after 20 weeks of gestation, complicated by at least one other endorgan dysfunction, such as proteinuria, renal dysfunction, pulmonary edema, and fetal growth restriction [1].PE infringes on 4.6% of pregnant women worldwide [2,3] and 2.3% of all pregnancies in China [4], and remains a major cause of maternal and perinatal morbidity and mortality [5].Regardless of the timing of onset, all PE patients are at risk of rapid deterioration into hemolysis, elevated liver enzymes and low platelet count (HELLP) syndrome or eclampsia [6][7][8], posing a huge burden to healthcare across the world.Except for these shortterm risks, the long-term impact, such as cardiovascular, renal, and cerebrovascular events, posed a threat to both the mothers and their offsprings [9][10][11][12].The etiology of PE remains largely unknown, but placental dysfunction is thought to be the core pathological change for PE.According to the two-stage model, placental dysfunction (stage 1), featured by oxidative stress, mitochondrial dysfunction, metabolism disorder, and apoptosis [13,14], is the initiating event for the subsequent maternal disease (stage 2) by releasing soluble toxic materials into the maternal circulation [15].During this process, multiple pathogenic processes, including angiogenic imbalance, immune cell activation, and endothelial cell dysfunction were involved [16][17][18], but the underlying molecular mechanism remains obscure.Currently, delivery of the fetus is the definite treatment of PE, however, preterm and the long-term risk of complications have brought impediments to PE management.Pharmaceutic therapies have been used in PE management.For example, aspirin is recommended for PE prevention in high-risk pregnant women [19].Other small molecules, such as statins and metformin, have been shown to reduce adverse outcomes and incidence of PE [20,21].However, there is still a great gap between the treatment options and clinical needs, especially in fetal outcomes.Therefore, a comprehensive understanding of pathogenesis and novel agent development will be of capital importance to PE treatment.Autophagy, an intracellular degradation system for damaged or dysfunctional cellular components, plays a capital role in maintaining cellular homeostasis [22].It is also the core molecular pathway that regulates embryo development during normal pregnancies [23].Intervention in normal autophagy processes will lead to a plethora pathologies, including cancers, metabolic disorders, and many others [24].Recently, increasing evidence suggests that autophagy is associated with PE [25][26][27], however, these studies provided conflicted results.For example, Akaishi et al. [28] reported a decreased level of p62, a substrate of autophagy, indicating activated autophagy in PE placentas.On the contrary, Nakashima et al. [29] observed the accumulation of p62 in placental samples from PE patients, suggesting inhibited autophagy in PE.Such a contrast may be result of a complex regulation network and the microenvironment that autophagy involves in.Therefore, a thorough understanding of correlations between autophagy and PE as well as its microenvironment may lead to finding novel targets for PE and the development of novel therapeutics.
Despite minimal treatment options for PE, development of novel therapies is ongoing, including pravastatin, Metformin, and aspirin.However, due to the high demand for safety profiles in pregnant women, safety equals to or even outweighs effect in drug development for PE.In China, medicinal herbs have been a major treatment option for thousands of years, including the pregnant women.Studies [30][31][32] accumulated that some active ingredients of herbs, such as berberine, baicalein, and luteolin, can exert a protective effect against PE.For example, berberine (BBR), a major component isolated from Coptis chinensis, is active on multiple organs in the human body.Studies showed that BBR can regulate metabolic disorders, reduce hypertension, and alleviate preeclampsia by regulating IL-2/IL-10 balances and inhibiting apoptosis [33,34].Baicalein (BCL), a flavonoid extracted from Scutellaria baicalensis Georgi, is known for its antitumor effect in various malignancies.It was reported that BCL reduced blood pressure by regulating inflammatory and vascular disease-related factors in hypertensive pregnant rats [35].Also, BCL was reported to protect liver and kidney function by inhibiting apoptosis in a PE model [36].As S. baicalensis Georgi is one of the herbs that were most frequently prescribed for and used by pregnant women [37], BCL seems to be a potential agent with a good safety profile for PE.In addition, luteolin was found to be a candidate for PE treatment by reducing anti-angiogenic sFlt-1 releasing and promoting vasorelaxation of uterine arteries [38,39].However, whether these active ingredients can target autophagy remains obscure.Therefore, this study aims to identify novel autophagy-related targets of PE and eventually find a safe therapeutic for PE by drug screening.First, a weighted gene co-expression network analysis (WGCNA) was applied to filter autophagic hub genes that related with PE.The differential expression levels and regulatory network of hub genes were also validated with two independent datasets and analyzed via specific toolkits, respectively.Associations between the hub genes and the infiltrated immune cells were analyzed to explore the role of autophagy in PE microenvironment.Additionally, as the potential therapeutics for PE, berberine, baicalein, and luteolin were screened by molecular docking in comparison to pravastatin, metformin, and aspirin.

Weighted gene co-expression network analysis (WGCNA) based on autophagic gene set
The dataset GSE75010 was downloaded from the Gene Expression Omnibus datasets (GEO, https:// www.ncbi.nlm.nih.gov/ geo/), which included 80 placental samples from patients with preeclampsia (PE) and 77 placenta samples from patients without preeclampsia (non-PE).The clinical characteristics of these pregnant women can be found in Table 1.Then, a list of 1334 autophagic genes (Supplementary file 1) was obtained from HADb: Human Autophagy Database (http:// www.autop hagy.lu/), AUTOPHAGY DATA-BASE (http:// autop hagy.info/), HAMdb (http:// hamdb.scbdd.com/), and GSEA (http:// www.gsea-msigdb.org/ gsea/ index.jsp).The R package WGCNA software was then used to construct a weighted co-expression network of autophagy-related genes that were differentially expressed in PE patients.During WGCNA analysis, we set the optimal soft threshold to be 8. Next, genes were clustered using average linkage and Pearson's correlation and a hierarchical clustering tree was constructed using dynamic hybrid cutting, with the minimum number of genes per module being 30.After determining gene modules by the dynamic shear method, the eigengenes of each module were calculated, and then clustered and merged the nearby modules into new modules, with a height equal to 0.25.Finally, we studied the module-trait relationships with Eigengene Networks.All data in this study were obtained from public databases and did not require the institutional review board's approval.

KEGG pathway and GO functional enrichment analysis
Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis [40], as well as Gene Ontology (GO) analysis, were performed to reveal the enriched pathways and biological functions of modules using the R package "clusterProfiler".Then, the ggplot2 package was applied to visualize the enrichment result.GO enrichment included biological process (BP), molecular function (MF) and cellular component (CC).

Hub genes screening and expression analysis
First, the gene significance (GS) and module membership (MM) were calculated.GS refers to the correlation between gene expression in the blue module and PE, while MM represents the correlation between the blue module and genes.Genes with |GS| > 0.5 and |MM| > 0.8 were selected into gene set 1.Then, a protein-protein interacting (PPI) network of the red module was constructed via the STRING database (https:// string-db.org/).The plugin cytohHubba of the Cytoscape software (https:// cytos cape.org/) was utilized to calculate the connectivity and selected the top 10 genes into gene set 2. Next, the common genes that overlapped between gene set 1 and gene set 2 were filtered to be hub genes.Finally, the Wilcoxon rank-sum test and Kruskal-Wallis test were used to analyze the expression of the hub gene in various groups and the correlations in between.

Validation of the hub genes with independent datasets
Two independent datasets, GSE24129 and GSE66273, were downloaded from the GEO database (https:// www.ncbi.nlm.nih.gov/ geo/).GSE24129 dataset contained 8 PE and 8 non-PE cases, while GSE66273 dataset included 6 PE and 5 non-PE cases.The difference in expression of the hub genes was validated with these two datasets via the Wilcoxon rank-sum test.The threshold was set as the adjusted P<0.05.

Potential molecular mechanisms of hub genes
The specific regulatory mechanisms of hub genes were explored with a multifactor regulatory network analysis.We investigated the interactions between hub genes and between hub genes and autophagic genes.Besides, the R package "clusterProfiler" was applied to analyze the signaling pathway that hub genes and autophagic genes involved.In addition, the networkanalyst database (https:// www.netwo rkana lyst.ca/) was used to predict miRNA-hub gene interactions, as well as transcription factors (TF) -hub genes interactions.Moreover, we studied the correlations between lncRNA and hub genes with the help of the starBase database (http:// starb ase.sysu.edu.cn/) and visualized it with Cytoscape software.

Infiltrated immune cells in PE patients
Infiltrating immune cells of the GSE75010 dataset were analyzed via the CIBERSORTX platform (https:// ciber sortx.stanf ord.edu/).LM22 was selected as signature matrix file and B-mode was enabled for batch correction to obtain the proportion of 22 immune cell types for each sample (Supplementary file 2).Then, the differences in the proportion of immune cells between the PE group and non-PE group were analyzed.

Drug screening via molecular docking for hub target proteins
The

WGCNA construction based on autophagic gene set
The flowchart of this study can be found in Fig. 1.Clinical parameters of the PE group and non-PE group are presented in Table 1.A total of 7 modules were obtained in WGCNA analysis.The correlations between the sample trait and module eigengene (ME) of each module were calculated.As shown in Fig. 2a, MEred correlated most to PE, which included 54 genes.GO enrichment analysis of these genes indicated that autophagy, regulation of autophagy, response to hypoxia, and cellular response to external stimulus were the major bio-function of these genes (Supplementary file 3) (Fig. 2b).KEGG enrichment analysis suggested that HIF-1 signaling pathway and TNF signaling pathway were mainly involved (Supplementary file 4) (Fig. 2c).

Hub genes screening and validation
To filter the hub genes of autophagy-related PE differentially expressed genes (DEG), we obtained 11 module hub genes (Fig. 3a) and 10 PPI hub genes (Fig. 3b), respectively.Then, a total of 3 genes, namely LEP, PKM, and HK2, concurred in both gene sets, which were considered the hub genes (Fig. 3c).It showed that the hub genes were differentially expressed in embryo samples of the GSE75010 dataset (Fig. 3d, e).Besides, two independent datasets also showed that these genes were differentially expressed between PE patients and the control group (Fig. 3f, g).

Correlations of hub genes and clinical parameters
The expression levels of PKM and LEP were higher in patients older than 35 years old than that younger than 35-year-old (p<0.05;p<0.05) (Fig. 4a).Besides, there were significant differences among the three BMI groups in the expression of hub genes (all p<0.05) (Fig. 4b).However, the differences in the expression of hub genes between the miscarriage group and the non-miscarriage group, and between the hypertensive group and non-hypertensive group were small, except that the expression of PKM was remarkably higher in the hypertensive group than that in the non-hypertensive group (p<0.05) (Fig. 4c, d).

Molecular docking results
The binding free energies of aspirin, metformin, pravastatin, and berberine to the target proteins were relatively low (>-7.5 kcal/mol) when compared to that of baicalein and luteolin (Table 2).As shown in Fig. 7, baicalein has the highest affinity to PKM (-8.5 kcal/mol), with 2 hydrogen bonds and 11 residues involved in hydrophobic contacts, while luteolin has the highest affinity to HK2 (-8.9 kcal/mol), with 1 hydrogen bond and 10 residues involved in hydrophobic contacts.

Discussion
In this study, we applied a WGCNA analysis based on public databases and performed GO and KEGG enrichment analysis for the most correlated module.It turned out that hypoxia and oxidative stress were the major bioprocesses in PE, which is consistent with our prior knowledge.After further screening via the String database, we found 3 autophagy-related hub genes overexpressed in the PE group, namely, PKM, LEP, and HK2.In two independent datasets, these hub genes were also upregulated in the PE group.We further analyzed the hub gene expression among different subgroups and found that PKM and LEP were remarkably upregulated in women over 35 years old than those younger than 35 years old.Besides, women with a higher BMI tended to present higher expression of hub genes, and women with hypertension had a higher expression of PKM compared to women without.These results indicated that aging and energy metabolism were associated with autophagy.PKM encodes a pyruvate kinase involved in glycolysis, which has four isoforms (M1, M2, L, and R) in mammalian cells.Among these isoforms, PKM1 is a constitutively active enzyme in normal adult cells [41], while PKM2 is conditionally active and predominantly expressed in fetus and cancer cells [42].Contrary to PKM1, the low glycolysis activity of PKM2 can limit glucose oxidation and allow cells to grow under hypoxia [43,44].This may explain the high expression of PKM (PKM2) in the placentas of PE group as a result of adaptation to the hypoxia environment.Therefore, activation of PKM2 might alter the low flux and improve energy depletion caused by hypoxia.LEP is located on chromosome 7, which encodes a hormone peptide leptin that plays a capital role in energy homeostasis.Except for white adipocytes, leptin can be also produced by other tissues, including endometrium, placenta, and umbilical cord, with a wide range of functions in regulation of angiogenesis, vascular function, blood pressure, immunity, and placenta development [45][46][47].
Studies reported that hypoxia could lead to placenta leptin gene expression and production [48], and abnormal placenta leptin level was associated with aberrant trophoblast proliferation or invasion, as well as endothelial dysfunction, hypertension and fetal growth restriction [49,50], which may lead to PE.The interacting relations between leptin and autophagy varies in a tissue-specific pattern [51].Although, the evidence in placenta is scarce, it was reported that autophagy was active in the placenta of PE women and leptin could modulate autophagy in various tissue types [52][53][54][55].HK2 is a protein-encoding gene that encodes hexokinase 2, an enzyme that catalyzes glucose into phosphorylation form.It was suggested that HK2 was pivotal to decidualization and non-coding RNA of HK2, including lncRNA and miRNA, could lead to dysfunction of HK2 and possibly cause PE [56].Therefore, we further investigated the regulation network among hub genes, lncRNA and miRNA to explore the underlying mechanisms of hub genes in the pathogenesis of PE and found that the signaling pathway involved were mainly autophagy and other metabolic pathways.These results indicated that autophagy might mediate PE via metabolism disorder.
Fig. 5 Multifactor interaction network of hub genes, autophagic genes, miRNA, lncRNA, and signaling pathways.Orange ovals referred to hub genes; green ovals referred to autophagic genes; blue ovals referred to lncRNA; yellow marks represented miRNA; red ovals stood for signaling pathways In addition to the placenta metabolic disorders, maternal factors may also be involved, such as maternal immune function.It is well established that immune-inflammation disorder at the maternal-fetus interface participates in PE development.For example, natural killer (NK) cells and Tregs in the decidua influence placentation; inhibition of NK response to MHC and deficient in Tregs may lead to impaired spiral artery remodeling [57,58].However, our studies showed an elevated proportion of activated NK cells and Tregs in PE group.This discrepancy might be result of mRNA-based analysis of this study, which only reflected quantities but lacked functional information of these immune cells.Macrophages obsess potential for phenotypic plasticity, which are featured by immunosuppressive M2 phenotype and inflammatory M1 phenotype.During normal pregnancy, M2 macrophages predominant in decidua to support fetus-maternal immune tolerance.Whilst in hypertensive pregnancy, aberrant macrophage activation may lead to insufficient trophoblast invasion and induce placenta apoptosis [59].In this study, the M1 macrophages in PE group was numerically higher than healthy women, but the difference was small; while M2 macrophages were remarkably lower in PE, indicating the imbalance of macrophages in PE placentas.Yang et al. [60] reported spontaneous miscarriage women displayed inhibited autophagy and reduced retention of decidua macrophages, which was mediated by a lipid metabolite.This study suggested that metabolism reprogramming can regulate immune balance via autophagy and subsequently result in poor outcomes in pregnant women.Multiple studies have mentioned the immune-regulatory roles of the hub genes.For example, PKM2 was found to promote M1 type activation of macrophages [61], which might explain the Fig. 6 The abundance of infiltrated immune cells in PE placenta samples and correlations between hub gene expression levels and immune cell proportions upregulated PKM and downregulated M2 macrophages in PE group.Besides, PKM2 could stimulate activation of NK and neutrophils [62,63], possibly leading to an inflammatory environment for PE pathogenesis.Nevertheless, these results were based on other tissue types.In the future, more studies should be performed on PE population to investigate the relationship between hub genes and immune cells in placenta.
Contemporary research has focused on inventive therapy for PE, including aspirin, statins, metformin, and other recombinant proteins.Aspirin is a nonsteroid anti-inflammatory agent that also reduces sFlt-1 overexpression to counteract PE [64].Similarly, statins and metformin were shown to reverse placental perfusion by sFlt-1 inhibition [65,66].In this study, however, these small molecules didn't display high affinities to autophagy-related hub genes according to the molecular docking results.Molecular docking is an essential tool for drug design with its ability to predict the conformation of small molecules and targets [67], and it has been applied in many studies to screen potential drugs [68,69].The docking results showed there was a low binding affinity of BBR to hub genes, indicating autophagy might not be involved in the working mechanisms of BBR when used in PE.In contrast, BCL and LTL displayed high affinities to hub genes with binding affinity being as high as -8.4 and -8.9, respectively, which means they may bind directly to the hub genes of PE and take effect.What's more, numerous research suggested BCL and LTL could modulate autophagy in various tissue types [70][71][72][73].For example, BCL can promote autophagy and reduce ROS in cardiac muscle [72], while luteolin can inhibit autophagy via mTOR signaling in lung tissues [70].Although the research rendered inconsistent results, which might be due to different tissues and diseases, at least they displayed a potential of BCL and LTC to target autophagy.Taken together, we speculate BCL and LTL might be the promising drugs that target autophagy in PE, but more solid proof should be provided in the future.
Limitations inevitably exist in this study.First, the results were based on public datasets which lack detailed clinical information about pregnant women.Second, the mode of action of hub genes and specific bioeffects of potential drugs after binding to hub genes remain elusive.More functional experiments and mechanism research should be performed to verify the results.In summary, this research has identified PKM, LEP, and HK2 to be promising biomarkers for preeclampsia, which might regulate the pathogenesis of preeclampsia via targeting autophagy, metabolism and immune microenvironment.Molecular docking simulation predicted baicalein and luteolin to be promising novel drugs for PE.Hopefully, this study might highlight the role of autophagy in PE and contribute to the drug development for PE.

Fig. 1
Fig. 1 Flowchart of this study

Fig. 2 Fig. 3
Fig.2WGCNA analysis for autophagy-related genes that are associated with PE. a WGCNA analysis for gene modules.b GO functional enrichment analysis of genes in MEred.c) KEGG pathway analysis of genes in Mered, which were prepared with the KEGG database

Fig. 4
Fig. 4 Differences in the expression of hub genes among subgroups.a Differences in the expression of hub genes between different age groups.b Differences in the expression of hub genes among BMI groups.c Differences in the expression of hub genes between different outcome groups.d Differences in the expression of hub genes among patients with or without hypertension

Table 1
Clinical characteristics of pregnant women a BP Blood pressure PE

group (n=80) Non-PE group (n=77) p-value
structures of the target protein encoded by hub genes were downloaded from the Protein Data Bank (PDB) (http:// www.rcsb.org) database, while the structures of the candidate molecules, berberine (BBR), baicalein (BCL), luteolin (LTL), pravastatin, Metformin, and aspirin, were searched from the PubChem database (https:// pubch em.ncbi.nlm.nih.gov).Then, the protein structure files were converted to *.pdbqt format and docked with drug molecules using AutoDock1.2.0 software by removing the original ligands and adding polar hydrogen atoms.AutoDock Vina software (ver.1.1.2,the Scripps Research Institute, U.S.) was applied to derive binding free energies to simulate binding affinities between receptors and ligands.Finally, the folding patterns and molecular interactions between target proteins and drug molecules, as well as the bonds within the ligand-receptor complexes were visualized by Ligplot software (ver.4.5.3,EBI, UK).

Table 2
Affinities of candidate molecules binding to target proteinsa Refers to number of residues involved in hydrophobic contacts